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1.  Introduction 

c 

"This  paper  consists  of  the  computation  of  several  hitting  time  and  hitting 
place  distributions  for  two-dimensional  Brownian  motion.  The  motivation  for  this 
study  is  two-fold:  first,  to  get  a  diffusion  model  for  the  firing  behavior  of 
a  simple  network  of  neurons,  and  second,  to  get  an  Interesting  two-dimensional  ver¬ 
sion  of  the  inverse  Gaussian  distribution.  '  /-• 

Fienberg  (1974)  has  reviewed  various  models  for  the  firing  of  single  neurons. 

A  classical  model  of  Gerstein  and  Mandelbrot  (1964)  says  that  if  the  electrical 
state  (or  potential)  of  the  neural  membrane  is  specified  by  a  single  number, 
which  moves  towards  or  away  from  the  firing  potential  as  the  neuron  receives  ex¬ 
citatory  or  inhibitory  input,  resp.,  then  the  time  to  firing  can  be  approximated 
by  the  first  hitting  time  of  a  certain  level  fox  a  Brownian  motion  with  drift. 

The  authors  showed  that  this  model  could  Be  used  to  provide  a  satisfactory  fit  to 
some  data  that  they  observed;  more  importantly,  they  showed  by  Monte  Carlo  methods 
that  neural  activity  in  the  presence  of  stimuli  could  also  be  well  duplicated  by 
a  modification  of  the  above  random  walk  model. 

Next,  the  review  by  Folks  and  Chhikara  (1978)  shows  that  the  inverse  Gaussian 
distribution  has  many  nice  statistical  properties  which,  to  a  large  extent  mirror 
those  of  the  Gaussian  distribution.  It  is  natural,  then,  to  ask  whether  there 
is  a  multivariate  inverse  Gaussian  whose  statistical  properties  are  similar  to  those 
of  the  multivariate  Gaussian. 

Several  proposals  for  a  bivariate  inverse  Gaussian  have  already  appeared  in 
the  literature.  Barndorff-Nielsen  and  Blaesild  (1983)  define  reproductive  ex¬ 
ponential  families  and  propose  a  bivariate  inverse  Gaussian  model;  they  claim  that 
their  generalization  has  nice  statistical  properties  (i.e.,  affords  tractable 
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estimation  and  analysis  of  variance),  but  their  proposal  does  not  have  inverse 
Gaussian  marginals.  Wasan  (1969,  1972)  proposes  several  bivariate  inverse  Gaussians 
but  does  not  develop  their  properties. 


A  Simple  Neural  Network 


Consider  the  three  neurons  of  Figure  1.  Neuron  A  sends  predominantly  exci¬ 


tatory  signals,  s,  to  B  and  C.  B  and  C  share  a  common  noise,  n,  and  they  also  have 


independent  sources  of  noise,  n^  and  n2>  resp.  If  the  electrical  states  of  B 
and  C  are  encoded  by  single  numbers,  X^(t)  and  X2(t) ,  resp.,  then  X^(t)  has 


three  components;  the  common  noise  n,  the  particular  noise  n^,  and  the  signal  s. 

2  2 

Let  the  noise  variances  be  0  for  n  and  for  n^;  then  it  is  easy  to  see  that 


corr(X1(t) ,X2(t)) 


id  +  -|)d  +  -|>r1/2 


which  is  a  function  of  the  noise  ratios.  Also,  we  may  allow  the  drifts  of 
X^(t)  and  Xj(t)  (due  to  the  signal,  s)  to  be  different  since  B  and  C  may  accept 
the  same  input  but  integrate  it  differently.  When  either  Brownian  reaches  the 
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firing  threshold,  the  appropriate  neuron  fires,  returns  to  its  resting  state, 
and  the  process  starts  afresh.  What  are  of  interest,  then,  are  the  firing  times 
or  the  first  hitting  times  for  the  Brownian  motions.  Alternatively,  this 
model  can  be  used  to  study  a  single  neuron:  if  we  postualte  that  the  neuron  has 
two  interacting  trigger  zones,  (Gerstein,  et  al.  (1964)),  then  the  components  of 
the  Brownian  motion  describe  the  electrical  state  of  each  zone.  Mathematically, 
the  two  problems  are  the  same. 


3.  Preliminaries 

We  start  with  a  correlated  driftless  Brownian  motion  X(t)  with  EX(t)  -0 
and  var  X(.t)  -  t$.  Here 


1  p' 

.P  1. 


and 


cos  8  sin  8 
sin  8  cos  8 


1/2 

where  p-sin(28),  |8|^_5  tt/4.  Thus  X(t)  *  t  Z(t)  where  Z(t)  is  a  standard 
Brownian  motion:  var  Z(t)  *  tl.  Also  define  the  two  stopping  times  T^-infit:  X^(t)  = 
*  inf{t:  Z(t)  Here  a^  >  0  without  loss  of  generality,  and  is  the  line 

(vclR  :v,#$  ei“ai^  an<*  ^ei»e2^  *s  standard  basis  for  1R  .  By  the  scale 
invariance  of  Brownian  motion,  we  can  take  a^  ■  1.  Finally,  by  elementary  methods, 
we  arrive  at  the  following  problem  (see  Figure  2):  start  a  Brownian  motion  at  x 
-  (xltx2)  ■  (rQcos  80  ,rQsin  eo),  and  study  the  stopping  times  and  places  associated 


with  and 
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Figure  2 


Here  is  the  first  hitting  tine  of  L^,  and  a*  -j  +  8^n  ^P*  Also,  let 
W  *  {(r,0):  r  >0,  0<  6<  a},  and  T*  *  a  be  the  first  hitting  tine  of  3W. 

Our  ain  is  to  get  the  joint  density  of  and  on  the  way  we  compute  other 

quantities  that  are  also  of  interest.  In  particular,  we  study  the  following 
quantities: 

a)  PX(t’>  t,  Z(t)  e  B),  B  c  W 

b)  Px(t'>  t) 

c)  PX(Z(t')€  A),  Ac  3W 

d)  PX(T*€  dt,  Z(t)  c  da)  ,  a  c  3W 

e)  PX(T1  €  ds,  t2c  dt) 

f)  the  above  quantities  in  the  presence  of  drift. 

Here  P  is  the  measure  associated  with  standard  Brownian  motion  starting  at 
x;  Ex  will  denote  the  corresponding  expectation.  Note  that  the  marginal  dis¬ 
tributions  of  T  are  easy: 

X  ^  A  A 

Px(t,j  €  dt)  -  — <)>(—  )  dt,  and  PX(t.  c  dt)  -  — -  )dt 

1  trt  /r  t/t  •'t 
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-1/2  -x  /? 

where  A  *  x^slncc -  X2<:osa  and  $(x)  =  (2tt)  e  , 


4.  Brownian  Motion  in  the  Wedge 

The  main  result  of  this  section  is  contained  in  (8).  Most  of  the  subse¬ 
quent  results  follow  from  it.  If  we  have  a  positive  bounded  continuous  function 
f  defined  on  W,  and  which  vanishes  on  3W,  then 

(1)  u(t,x)  *  EXf  (Z  (j’a  t))  *  [  f  Cy)PX(f>  t,Z  (t)  €  dy) 

JW 

satisfies  the  heat  equation  with  boundary  and  initial  conditions 

(2)  ut  Au  in  W;  u(0,x)  *  f(x),  x  e  W;  u(t,z)*0,  z  e  3W. 

We  can  solve  C2)  when  a  =  Tr/m,  m  =  l,2,,..  by  the  method  of  images.  That  is, 
if  we  let  Tq  - 1 ,  F^  be  the  matrix  representing  the  reflection  across  the  line 
y  *xtan  ,  and  T^  -  F ^  o  T^_^,  and  let  f(y)  -  (-1)  f(T^y)  for  ye  T^CW),  we  have 
the  initial  value  problem 

1  ~  ?  ~  2 

(3)  ut«2^u  in  IR  ;  u(0,x)»  f(x),  x  e  1R  . 

The  solution  is  ' 


r  k  X-V  dv 

f(y)  l  <-i>V r->  -f 
JW  k-0  /t 


u(t,x)  -  E  f (Z(t)) 


where  4>2(x)  *  (2tt)  exp(-x'x/2).  It  is  easy  to  see  that 

.  2m- 1  .  x-T.  y 

(4)  PX(i’>t,Z(t)  e  dy)  -  i  [  (-1)V( - — )dy. 


While  (4)  is  appropriate  only  for  the  special  angles  ©Q  -  w/m,  the  following 
argunent  gives  us  the  result  in  general.  To  facilitate  this,  we  use  polar  co- 
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ordinates:  let  x  -  (rQcos0o,  rgSinOg),  y  =  (rcos0,  rsin©)  to  get 

2  2 

r  +ta  rr 

1 - k  ~~t~ 

(5)  PX(T'>  t,  Z  (t)  e  dy)  «  -r—  e  I  (-1)  e 

Z1Tt  k=0 


where  0fc  is  the  argument  of  Tky.  Rote  that  (Magnus,  et  al.  (1966)) 


iyz  -  2  y  T  (Y)I  (z) 
'■  n  n 


where  T  is  the  n^  Chebyshev  polynomial  and  I  is  the  modified  Bessel  function 
n  n 

of  order  n.  Recalling  that  T  Ccos0)  -  cos(n0),  and  -using  the  fact  that 


y  (-1;  cos  n(0-0,  )  •  2m  sin(n8)sin(nBg) 
k*=0 


if  m  divides  n  and  zero  otherwise,  we  get 


2 ,  2 
r  +r 

. - — ^  oo  mr  0-  rr 

(8)  PX(.T*  >  t,  Z(t)  €  dy)  — r|  e  1  l  sin  sin  — —  Intr(~T'’)drde 

”*°  T 


whenever  a  *  ir/m.  But  formula  (8)  is  valid  for  all  a,  and  it  is  easy  to  see  that 
it  solves  problem  (2). 

Expression  (8)  has  also  been  essentially  derived  by  Sommerfeld  (1894)  by  an 
extension  of  the  method  of  images.  See  also  Carslaw  and  Jaeger  (1959)  and 
Buckholtz  and  Wasan  (1979). 

We  next  compute  PX(l*>t).  If  we  integrate  out  0  and  r  in  (8)  and  use  the 
following  identities: 


2  I\(*)  - 


r  e'6t  iv(at)dt  •  i\/f  e"a  /8B 

i  n 
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(Magnus,  et  al.  (1966))  we  get 

2 


3r 


(10)  Pxct*>t) 


ifa; 


_0_ 

4t 


'Tit 


nil 
n  odd 


1  .  OTeo 

—  cin - 

n  a 


Xv+l(4t)J 


where  v  ■  off/ a. 

When  the  wedge  angle  is  special-  a-  ff/m-  we  have  the  alternate  expression 

2m- 1  ...  . 

(11)  PX(f>t)  «  2  £  (-l)k+1F(^>. 

k-0  “ 

where 

fu  rn  rn 

(12)  F(u)  -  cosC0-eft)R(-  -S’  cos(B-ert))d6. 

Jo  /t  °  Jt  ° 


Here,  R  is  Mills'  ratio:  R(x)  «  $(-x)/$(x),  where  4>  and  <j>  are  the  normal  dis¬ 
tribution  and  density  functions,  reap.  We  omit  the  details  of  this  computation. 

The  quantity  PX(t  >  t)  was  also  studied  by  Spitzer  (1958) ,  who  computed  its 
transform.  Checking  the  asymptotics  of  modified  Bessel  functions,  it  is  easy  to 
see  that  EXCflP»  f  T^-1PX(l*  >  t)dt  <  °°.  if  arid  only  if  crf3<  tt/2,  independent  of 

Jo 

x. 


The  distribution  of  ZCt1)  is  also  of  interest.  Now  u(x;A)  «  P 
isfies  Laplace's  equation  Au  -  0  with  boundary  condition  u(x,A)  « 
The  Green's  function  for  the  wedge  is  easily  computed,  and  we  have 


(13) 


-  l 

'  &)a  "  sin 

Px(Z(T')eda)  -  ±  - 


ire. 


sin2  +  (1  +  cos  -~)2  r° 


x(Ztf)  e  A)  sat- 
l{x  e  A}. 


where  we  use  the  plus  sign  for  T2  <  and  the  minus  sign  otherwise.  Using 
elementary  estimates,  it  is  easy  to  see  that  ExZi  (j*f  exists  iff  a8<ir,  again 

*  # 

Expression  (10)  corrects  a  mistake  in  Wasan  and  Buckholtz  (1979). 


8 


independent  of  the  initial  position,  x. 


5.  Joint  distribution  of 

Using  an  argument  similar  to  that  of  Daniels  (1982)*,  it  can  be  shown  that 
if  PXCt  >  t,  Z(t)  € dy)  «  f(t,x,y)dy,  then 

<14>  PXCt,€dt,Z(t)  eda)  -  fCt,xfy)|  Jdadt 

g 

where  "3n  denotes  the  derivative  in  the  inward  normal  direction.  See  Figure  3. 


.x 


it-  n 

V 


/ 


S  An 

Y— X- 


Figure  3 


When  6  ^  0(cH,  ~  ~  so  from  C8)  we  get 


(15)  PX(T  edt.ZCr’Je  da) 


TT 

2. 

a  ta 


2  2 

a4To 

2t 


00  mre. 

J  i  n  sin - — 

n-0  n  a 


ar 

1  -(— ) 
nTT  t 

a 


where  6  is  1  if  0  ■ 
n 

PX(T'€  at,X(T')€  da)  - 

the  line  y  ■  x  tan  ~. 
able: 


0  and  (-1)  if  8  ■  a.  It  is  clear  by  symmetry  that 

X 

P  (1  e dt,Z(T’)€  da’)  where  x  is  the  reflection  of  x  across 
And  for  the  special  angle  a  ■  Tr/m,  a  simpler,  formula  is  avail** 


* 

I  thank  Professor  D.  Siegmund  for  this  reference. 
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2 0-1  k  ,  ,x"Tuy, 


(16)  Px(t'€  dt,ZCT')€  da)-  (-1)K  4>2(— — ^-) 


/r 


with  Px(T,e  dt ,  ZCrOeda')  computed  by  symmetry. 


y,=o 


<x’Tke2> 


Finally,  we  can  compute  Px(t  <z  ds ,  x  e  dt) ,  the  joint  density  of  (x  ,x  ) 

12  1-2 

By  the  strong  Markov  property  we  have  for  s  <  t. 


(17) 


P  (x,  e  ds,x9  e  dt) 


Px(x1  e  ds,  ZCt')c  da,  x2  e  dt) 


3W 


I; 


PXCx1  =  t*  e  ds  ,  Z (i  ‘  )e  da)  P®  S±n  a(x2  c  dt-s)  . 


But  the  first  term  of  the  integrand  is  given  by  (15)  and  the  second  term  is  just 


(18) 


,a  sin  a 


(x2  c  dt-s)  -  3  ----- 


since  x^  is  just  a  one-dimensional  inverse  Gaussian.  After  some  computation 
and  (9)  we  get 


(19)  P  (x^  e  ds,x2  e dt)  - 


tt  sina 


mre. 


.  0  2(t-s)+(6-s  cos  2a)  r-r  .  iviri-l  0_  r 

TT-.l+frU  2al  1  i. t_U  "“"V1./' 


rQ  Ot-s)/s 


2o2/r(t-s>C - '  -  2a  a's)+(;t-s  cos  2“>  n-0 

t— 8  cos  a 


nTT  I~2  Ct-s)+2  (jt-scosl 
2a 


Finally,  using  the  fact  that  (see  Figure  3) 

PX(X1  c  ds,x2  e  dt)  -  Px(xx  e  dt,x2  c  ds) ,  s  <  t 

the  joint  density  of  (x^,x2)  is  determined  for  s  4  t.  In  (19)  we  can  let  t-*-s  and 
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use  the  fact  that  as  z -*■(),  I^(z)  ~  (-|)V/r(v+l)  to  get 


(20)  PX(T e  ds,x2  e  dt)  + 


4ns 


r0sin260 


exp(-  -£) 


if  0  <  o<- 
2 

if  --  <  a  <  tt 
2 

if  a-  it/2 


Thus  the  joint  density  of  (t^t  )  is  discontinuous  on  the  line  s  «=  t  only  when  the 
original  Brownian  motion  X(t)  is  positively  correlated.  Of  course,  we  could  have 

started  with  (16)  to  get  a  simpler  expression  for  the  joint  density  of  (T1»T2) 
for  the  special  angles  a  *  ir/mj  we  omit  the  straightforward  calculation. 


6.  Brownian  Motion  with.  Drift 

Of  course,  the  case  of  Brownian -motion  with  drift  is  of  more  interest.  The 
analysis,  however,  is  considerably  more  complicated  and  so  it  is  not  always  possi¬ 
ble  to  evaluate  the  integrals  that  arise.  In  this  section,  we  extend  the  results 
of  the  previous  sections  to  Brownian  motion  with  drift. 

If,  in  section  3,  the  correlated  process  X(t)  has  drift  6“  where 

6^>0,  it  is  easy  to  see  that  for  the  corresponding  process  Z(t),  in  the  wedge  W, 

2  1/2 

we  have  drift  p  *  (p1,p2)'  *  ^2”®!'  _02  (1-P  ).  4  Let 


(21) 


f Ct;a,b) 


exp  C- 


_  Ct-a)4 


2a 


be  the  inverse  Gaussian  density  in  its  usual  form  { 6  ,p.  263J.  Then  it  is  easy  to 

1*1 

see  that  ^  has  density  fCt;  - 


x.sina-x  cosa|  2 

— — — ■  ,  (x^sinot-x^osa)  ;  and  t2  has  density 


fCt;  ~=Z> 


Vl-P2 
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Let  P^  be  the  measure  associated  with  uncorrelated  Brownian  motion  starting 
at  x  and  with  drift  y; 


C22)  P^  (Z  (t^)eA^, •  •  • » Z  Ctn)eAn)  =  Pq  CZ  (t^)+y  t^  ,Z  (t^)+y  t^e  A^) 

for  all  n,  for  all  t,<...<t  ,  and  for  all  Borel  sets  A..  Our  basic  tool  is  the 

j.  n  i 

exponential  (likelihood  ratio)  martingale 


(23) 


dP 


dP 


=  exp  (y  ’  (Z(t)-x)-t  |y  |  /2) 


on  Ft»  the  sigma  field  generated  by  (Z(r)  :  s^t).  Thus,  we  have  that 


2 

(24)  vCt.x.y)  -  P*(t'  >t,  Z(t)€dy)  -  ey'  Cy-ac)" I  t/2P*(ib  t,  Z(t)e  dy) 


is  a  solution  to  the  diffusion  equation  with  convection  or  drift: 

(25)  vt  -  ■|^v+y,^v;  v(0,x,y)  -  6x_y;  v(t,a,y)  =  0,  a  e  3W 

where  6  is  the  Dirac  delta.  Of  course,  in  (24),  Pq(t'  >  t,  Z(t)e  dy)  is  given  by 
(8).  The  expressions  for  P^Cr'  >  t)  and  P*  ( Z(t ' ) e A)  do  not  seem  to  be  convenient 
as  they  were  for  the  driftless  case  ((10)  and  (13)  in  section  4),  but  the  joint 
density  of  T*  and  Z(t')  is  available.  In  fact,  Daniels’  argument  gives 
(see  Figure  3) 


(25)  P*(T 'cdt ,Z(t* )cda)  -e“X,y_t'V'  /2+VJ'a  P^x’edt.zCr'Ma) 

P*(.T-cdt,z(i')£da')  -  e-=‘Vt|v|2/2+.-()i1coS»Vlna)p*(T,td£t2(T,)e<i4,)> 


where  we  use  (15)  for  Pgfx’edtjZfr' )  cda). 
Finally,  we  have  for  s  <  t 
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?*(T1  e  ds,T2  e  dt)  -  j  P^(Tieds,2  (T^ 


,  Sna  sin  a.  x 

£da)P  (t2  £  dt-s) . 


Note  that  P^  f in  a(T2  £ dt~s)  is  just  the  inverse  Gaussian  density  and 

P*(t^  eds,  z (t^)  £  da)  is  given  by  (21).  The  integral  involved,  however,  does  not 

seem  to  be  tractable. 


7.  Concluding  Remarks 

Clearly,  the  y=  0  case  is  much  easier  than  the  y^O  case;  however,  the  physical 

motivation  demands  y ^ 0.  For  higher  dimensional  problems,  similar  methods  can 

be  used  to  get  the  joint  distribution  of  T«  ),  where  T.  «=inf{t>0:  X.(t)=a.}, 

i  p  i  i  i 

a^  >  0  ,  and  X(t)  ~  N(0,t$)  is  a  driftless  correlated  Brownian  motion.  The 
geometry  in  ]RP  is  quite  complicated  though;  for  certain  patterened  X  (e.g., 

X^  = p) ,  the  transformation  to  independence  can  be  done  in  closed  form,  and  the 
joint  distribution  of  T  is  available.  The  same  problem  for  general  $  and  for 
y jf 0  requires  further  investigation,  and  will  be  the  subject  of  a  subsequent 
paper. 
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